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Abstract 


A numerical method, which is simpler and more efficient than others currently in use, is proposed for the 
computation of the full viscous flow over an aerobrake body in hypersonic stream at high altitude. It treats the 
shock layer surrounding the blunt forebody and the near wake behind the base simultaneously by formulating 
the Navier-Stokes equations in conformal and azimuthal-angle coordinates. The computational domain is 
confined by the body wall, outflow surface and the bow shock, which is adjusted along the coordinate normal to 
the wall in the course of iterations. Because of the optimal grid and a well-developed alternating direction 
implicit factorization technique for the governing equations, reasonably accurate results can be obtained with 
a 28 x 36 x 7 grid and 400 time-marching iterations. Excellent agreement of shock location is found between 
the present result and the schlieren photograph. Details of the base flow and shear-layer impingement on the 
cylindrical aft body are presented for an adiabatic- wall case. 

L Introduction 


Recently, urgent reasons have arisen for finding the complete fiowfield solution about a finite-size blunt 
body flying hypersonically at high altitude. Since ground facilities cannot adequately simulate the fiowfield 
and enable determining the aerodynamic and heating loads needed for designing a spacebound vehicle 
traveling between low Earth and geosynchronous Earth orbits, analyses will depend on the numerical solution 
of the model equations. Although computational capability is available to treat three-dimensional (3D) flow 
over aerodynamic configurations at transonic speeds, straightforward applications of the available method- 
ology to hypersonic flow would require prohibitive amounts of computer time and storage. The favorable 
approach is to create a grid such that the bow shock and the body become coordinate lines, and the gradients 
along the body can be properly resolved. Along these directions, a simpler and faster method is proposed to 
facilitate the grid adaption to the wall and to accelerate the convergence rate of the numerical algorithm. With 
this method, a reasonably accurate viscous flow solution can be obtained from the full Navier-Stokes (NS) 
equations on a coarse grid after a few hundred iterations. Of more importance, its performance is equally 
impressive for flow at high angles of incidence and for arbitrary body shape. 

In the hypersonic regime of interest, the turbulence model and geometry details may be considered 
relatively unimportant compared to the dissociative nature of the air due to high temperature, yet numerical 
calculations of a finite-body flow are still limited in scope by the computer time and memory and by the 
numerical algorithms for equations and grid. Nearly all computations have been based on the assumption of 
laminar flow without finite-rate chemical reactions, and most are for two-dimensional (2D) flows. Lombard et 
al.(l>2) US ed a Navier-Stokes formulation on generalized coordinates and simple algebraic grid-clustering 
equation to study the fiowfield of the Galileo probe and of an aerobraking orbital transfer vehicle. The equa- 
tions are cast in curvilinear coordinates, but the vector variables are defined in Cartesian coordinates. In Ref. 

1, an alternating direction implicit (ADI) factorization technique was used, whereas a Gauss-Seidel type 
relaxation technique was used in Ref. 2 with nearly an order of magnitude reduction in the number of itera- 
tions. Gnoffo(3) has introduced an adaptive-grid algorithm to attack similar problems, but used a finite volume 
and the MacCormack explicit scheme. The grid distribution as determined by weight parameters had en- 
hanced the solution accuracy without adding more grid points to the computational domain. The techniques 
used in both Ref. 2 and Ref. 3 have undergone continual modifications to improve the accuracy and to reduce 
the computation time but are restricted to axisymmetric flow currently. 

An alternate method was proposed by the author to solve the 3D equations in curvilinear-orthogonal 
coordinates by a proven ADI factorization technique. (4) The full NS equations are cast in the conformal 
coordinates defining the planes at constant azimuthal angle. The transformation of a family of ellipsoid-cones 
resembling the Apollo command module and the Viking aeroshell probe onto a near half-circle is made possible 
through one application of the Karman-Trefftz equation as suggested by Moretti.(S) Thus, the difficult problem 
of initializing the complete flow over an arbitary body becomes a much simpler problem for a sphere. Further- 
more, the analytical mapping of coordinates ensures a smaller size of discretization error in a largely deformed 
grid in comparison to the grid generated by differential equations. The pros and cons of the classical mapping 
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versus other grid-generation techniques as well as the idea of combining two-dimensional orthogonal grids on 
slices through the three-dimensional field have been discussed by Ives.(6) 

In the present paper, a typical aerobrake configuration consisting of an ellipsoid, a flat disk, and a 
cylinder is used to demonstrate the capability of the new method. Details of the various coordinate transfor- 
mations and mappings are described in Sec. 2 of this paper. A derivation of the full NS equations in computa- 
tional coordinates appears in Sec. 3, followed by a description of computational boundaries in Sec. 4, the 
definition of initial conditions in Sec. 5, and a discussion of the numerical algorithm in Sec. 6. Results are 
presented for the aerobrake at an angle of attack of 1 8°. 

2. Coordinate Transformation 

The objective is to map the physical domain enclosed by the bow shock, the body, and the plane Z = 
constant onto a parallelepiped such that the boundaries become part of the computational coordinates (Fig. 1). 
The basic fluid-dynamics equations given in orthogonal curvilinear coordinates are used because of their versa- 
tility and simplicity to handle the dominant features of the hypersonic flow and a wide variety of body configu- 
rations. The equations are expressed in terms of the orthogonal coordinates £,q,£ and their scale factors 
h v h 2 ,h 3 in the direction of increasing $,q,C parallel to unit vectors$$&. The change in the position vector r can 
then be written as 


d?= / ih l d $ +fy 2 dt\ +^h 3 d< (1) 

in which both ($#) and (h v hM are functions of (5,q,Q. The physical space used in this study is defined by 
the spherical coordinates r,9,(p shown in Fig. 1. They are related to Cartesian coordinates X,Y,Z by 

X = qsin^cosC Y = q sinl^sin 4 Z = — qcos£ (2) 

where $ = n - 0, q = r, and £, = <t>. 

Then, conformal mappings are introduced to transform a given geometry in the Y = 0 plane to a near 
circle in the Y = 0 plane. An axisymmetric body can be easily obtained by rotating the circle image around the 
Z-axis such that each t, = constant plane has the identical image. The mapping is represented by 

l = Id q), q = q (i q), l = £, (3) 

The inverse of Eq. (3) is applied to find the physical image from a given shape in the mapped space. The first 
use of conformal mapping for numerically analyzing axisymmetric flow over a blunt body was by Van Dyke(^) 
in his well-known inverse method, by which both shock and body are conformed to coordinate lines with a 
single analytical function. 

In a direct manner to solve for the flowfield, the body shape is given in physical space, then transformed 
onto a smooth contour image. A general mapping function, known as the Karman-Trefftz equation, is applied 
repetitively to remove the corners and edges of complicated geometry.(S) It is exemplified by 

x~ h _ ( l - 1 W (4) 

X+h* 

where h and 8 are parameters, x and x refer to the complex function Jt = ^ + tq and £ = S + iq. The conjugate 
function of h is h*. The Joukowski equation is a simplification of Eq. (4) by choosing 8 = 0.5 and h = h* = 

R(h). 

The Cauchy-Riemann conditions stated in polar coordinates are 





n*. 


= - q^n = i 
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where R and I are the real and the imaginary part of the functional gxl X, which is composed of 


„ — _ — nJ® „ — „ jl v — r 
g - — - be , X - ne , X-ne 

dx 

Note the conformal property implies that ^ ^ + w = 0; hence, Rigtf j) — G and I(gx^x) = 0- The following 
relations will be used frequently. 


^ = Gn/n, n n = G, ^ n = o, n ? = o (£ 

The coefficient of magnification G and the angle of rotation co vary from point to point. 

The change in the position vector r now becomes 

5r=^ k /i 1 d £ d q + ^/i 3 d\ ^ 

The scalar factors in the mapped space are = q /G, h 2 = 1/G, and = q sin The next step is to 

make use of the shear transformation to bring both the outer and the inner boundaries of the domain of interest 
depicted in the mapped space into becoming a part of the coordinate lines in x,y . This transformation, shown in 
Fig. 2c, is accomplished by 


x = I/U, y=(S-q)/(S-B),* = h 2n 

where B = S = are the body and shockJocations_in the mapped space. More algebraic 

transformation may be used to cluster points toward q = B or $ = t and to result in a nonuniform grid in Fig. 
2c, but uniformly spaced grid in the computational space, Fig. 2d. The transformation equations used are as 
follows: 


x = 


c -f arcsinh f xl x 


— 1 ^ si 


sinh (coc)/co 


/ x 


c— In 


1 4* (exp(co) — 1) x Jx / 1 + (exp( — co) — 1 )x J x 


/ 2co 


y = In 


(p + y)/(p - y)j /in [(p + i)/(p - i) 


where (3 and co are used to control the degree of clustering. 

The transformations between the computational and the physical space are summarized schematically 
in Fig. 2. The mapped coordinates are given by_Eq. (4); the velocity components at each point are decomposed 
on a local polar frame which is rotated by co = £ — The computational coordinates are determined by Eqs. 

(7) and (8); velocity components are identical to those in the mapped space. 
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3, Navier-Stokes Equations in Computational Coordinates 


The conservative law of a single-component, compressible, viscous flow is commonly written in terms of 
orthogonal curvilinear coordinates 5,q,5 or 5,q,C To recast the governing equations in terms of computational 
coordinates, repetitive applications of chain rule are needed to result in the following form. 


u+f+g + h + r= 0 
t *x °y z 


( 9 ) 


where u = h 1 h u 


f= h 2 h Z X l f ’ S = h \ h #l f + h Z h i y Z 8 + h l h #l h + *1*2*3 y t U 


h = h^h^-h. 


0 


MW( 3 ) - h 3 (h 2 )-g( 3) -h 2 (h 3 )jh( 4) 

MW( 3 > - *i(*3*;*( 4 ) 

^(*3 

0 
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The scalar factors are defined in Eq. (6), and the metric coefficients of the nonorthogonal coordinate transforms 
are analytic. Convective and diffusive fluxes in the NS equations have standard notation. The velocity com- 
ponents u, u are equal to 


(“) = ( cosco - sin co\/«\ where co = 5-5 
\ v ) \ sinco cos co / \ u ) 

Other variables as well as the stress tensor and heat-flux vector can be found from standard texts. 

A set of simplified equations deduced from Eq. (9) by eliminating h 3 is used on x = 0 and n, as follows. 

U t + fx + 8y + r=0 

where u = h l h 2 u, f- (2 -u*/u) h^f + h z ~ uu*, g = h gr-f + h pr-g + h { hg t u 
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r = 


8 ( 1 ) 

8(2) + (hj-fo) - (h 2 ) X - i(3) 

8(3) + (h 2 )lf(3) - W(4)/A 3 - (A,)-^2) 
0 

8(5) 


The derivatives of (ft 3 )| and (/i 3 )“ are obtained by making use of the chain rule; hence, 

( ft 3^= = '■i^^nd M" h 2 sin{ ' 

Other derivatives of h t and h 2 , however, must be numerically determined using the following equations. 

(). = ,.( )x + y_ ()y ()- = y-() y 

4. Treatment of Boundaries 


The outer boundary consists of three segments, denoted by 1, 2, and 3 in Fig. 3. Segment 1, ranging from 


0 ^ ^ ^ £ max , is the bow shock, which varies in position according to the equation 


6s k+1 =bs k 




(ID 


where S * is estimated at the end of each iteration k by matching the Rankine-Hugoniot (RH^results to the 
flowfield results. When the computation converges, the shock positionjvilljiot change since S t =0. Segment 
2 represents a cone with a half-angle determined by the shock slope at £ = £ ^ The cone terminates at the 
plane Z = Z^, a specified parameter to control the extent of computation region. The radial location of seg- 
ments 2 and 3 is given by simple equations. The conic shock wave has no speed, so the RH equations alone are 
used to define the flow variables. The outflow surface, segment 3, has inviscid outward supersonic flow for its 
major portion; two-point extrapolation formula is considered adequate for prescribing flow variables therein. 

Let S* be the shock location along q at {* = n - 0 majc ; the cone is governed by the following equation. 

S = S* (sin e* - cos 0*tan( 0 + 0* - n/2)) / (sin 0 - cos 0 tan ( 0 + 0* - n/2^ 

The outflow surface is given by 

S = Z^sec 0 where©' = arctan ^ S * + Z^tan 
The corner of segments 2 and 3 is given by 


0 cr = arctan^ S * sin 0 * — (z S * cos 0 tan ^ 0 + 0 * — n/2 j j/Z ^ 

The velocity vector on the boundary first is computed in (£,q,0 components from Eq. (9), then is adjusted to 
satisfy the RH equations at the bow and conic shocks, and to meet the vanishing normal or no-slip conditions at 
the wall. In either case, an intrinsic coordinate system is needed to relate the local normal, tangent, and 
binormal vectors to the physical orientation vectors U$£). Let the body and the outer surface be 6 = q - 

B($,£) ands = q - S(£,0, respectively; the corresponding normal vector £ is (-B^/B, 1, -B^/(B sin £,)) and (S^/S, -1, 
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SMS sin $)), which are meant to point inward the computational domain. Other components of the intrinsic 
frame are obtained from£= n X K , and o = IX n. The procedures of shock fitting and wall velocity adjust- 
ment are identical to those described in Ref. 4, with the exception that u,v are obtained from u,v from 

(“W ®° S0) sin co V u\ where(0= ^ 

\U/ \ -sm 0) COSO) A y / 

After the boundary conditions have been imposed, an inverse transformation is used to get u y v. 

Velocity conditions may be imposed on the boundaries using the mapped coordinates. The unit normals 

to the shocks and the body are, respectively, equal to - BJ h {f 1 / /* 2 , — B {Sjfh l9 — l/h^S 
The zero body normal velocity is automatically satisfied in the physical space. Upstream, at the shock surface, 
the free-stream velocity components are cast in terms of the mapped coordinates by inverse transformation. 

5. Initial Conditions 


Two sets of flow variables are specified along the outer and inner boundaries. An initial estimate of 
shock has two parts. 

S = B + 8 (l + (If) for l S l cl and S = B + 8 (l + d\) for \ > l c[ 

where.cZ and e are input parameters for specifying shock shape; 8 is the standoff distance. Using the RH equa- 
tions, shock properties are obtained consistent with the slope. On the body, the Newtonian pressure and the 
isentropic equation are used to prescribe variables for £; < and base flow approximations wall density p w = 
0.25P*; energy e = CT^ where C y is specific heat at constant volume and T is temperature; velocity compo- 
nents^ = v = w = 0 are applied to the rest of the body, T, > is the body point where the body normal first 

turns into the free-stream vector starting from the stagnation point. The outflow surface has properties 
varying from the shock values at the conjuncture point between segments 2 and 3 to the free-stream condition 
at fj = n. Once the two lines are known, a simple interpolation is used to estimate the initial values in the com- 
putational space. The same procedure is valid if physical coordinates are used in the computation. 

6. Numerical Algorithm 

The computation space ( x,y,z ) is discretized into a three-dimensional network of grid points, which are 
defined by 


x m — (m - 2) Ax, y n — (n - 2) Ay, z e = (£ - 2)Az 

with Ax = l/(mc - 2), Ay = l/(nc- 2), A z = l/(Zc- 2), and m = 1, 2, ... me, n = 2, 3, ... nc, andf = 1, ... Zcp, where 
lep = lc + 1. The subscripts m and t cover a range of grid points beyond the computation domain in order to 
incorporate the boundary conditions. The relations indicate that any point in the computation space can be 
addressed by the intersection of three independent families of coordinate lines, but the metrics and the depen- 
dent variables still require three-dimensional arrays such as p n ^ = p(x m ,y fl ,zp. 

The resultant finite-difference counterpart of Eq. (9) is 

8 > + V + 6 / + V 1 + %,.< = 0 <l: 

where 8 1 + denotes the two-point operator, and 8 X , 8 y, and 8 2 are centered spatial operators except along n — nc 
at time level k + 1. 

The spatial approximations are solved simultaneously with the unknown vector u n m ( to allow greater 
changes of 8 t + u and, hence, to result in a faster rate of convergence from initial approximations to the final 
solution. A factorization technique due to Beam and Warming(8) which emphasizes the noniterative nature of 
the well-established ADI method and solves for the incremental vector of the conservative variable u has been 
adopted and modified. In this version, because the primitive variables v = (p,u,v,w,e) T are used, the Jacobian 
matrices have simpler expression and involve less arithmetic. The solution procedure comprises five steps. 



(j + S A - 8 JP + «i,ny w = tf - : 1 - <A«.) 

( + 8 B - 8 (£ + xrf,/) ) Au** , = tAu* . 

\ y yy / / rc,m,c n,m,£ 

(iI+SC-5 (F + xd,7)jAw fe + i = xAi;** . 

\ z zz I J n,m t £ rc,m,c 


(13) 


where x is the inverse of the local A t multiplied by a constant Courant number, and dj and <f, are the implicit 
and explicit damping coefficients. Each of the three equations in the middle of Eq. (13) represents a tridiagonal 
system of linear equations with 5 by 5 matrix coefficients. 

It is seen in Eq. (13) that the first step is to find the residual explicitly, then to smooth or to filter out the 
k 

localized errors in A u ^ by subsequent implicit steps performed along each coordinate. Obviously, the impli- 
cit calculations are noi’ needed for numerical stability consideration alone if A t n m t meet the Courant- 
Friedrich-Lewy (CFL) criterion, or any step can be deleted if (A^) min > CN(At^ m ^) minf where i— n, m, or t 
and CN denotes the Courant number. Nevertheless, the primary function of implicit solutions seems to be 
more in the acceleration of the overall rate of convergence than in maintaining solution stability. Numerical 
experiments have suggested that implicit solution is even preferable (without concern in computation cost) to 
the explicit solution because the effect of boundaries can be transmitted across the entire line. 

The Jacobian matrices of the convective and diffusive fluxes with respect to v are given as follows. 


A u = P Ay 

A f — PA Ay + PD Av x 
A g = PB Ay + PE Av y 
A h = PC Ay -1- PF Ay 

z 


where P = ^^2^3 


1 0 
a p 
v 0 
w 0 


0 

0 

P 

0 


0 

0 

0 

p 


q p u py pty p 


and q = (a 2 4- y 2 + u> 2 )/ 2 


The components of the Jacobian matrices are as follows. 


^ll — ux ^h l + wx^/hy - pxJh v A 13 P x^/h^, A 2l gex-/^ph^j, A^ A^, A lg (y *) x ^i 

A 33 ~ A ll’ A 41 = ^®y(p*i)> A 44 = A U’ ^46 = ( Y “ 1 ) A 52 ~ A 54 = px {/( ph 3 )’ ^55 " A n 

B 11 = °^| + V yJ h 2 + Wy i! k V B 12 = P yjj h V B 13 = P ^ ft 2' B U = P^ /A 3> B 21 = P*l) 

B 22 = B U' B 25= (v B 31 =^(P h 2 )' B 33 = B U> B 35 = (* " 1 ) ^*2’ B 41 = 

B 44 = B ll’ B 45 = ( Y “ 1 ) V^3’ B 52 “ Py ^( P ^l)’ B 53 = Py ^( P ^2)’ B 54 = W^(pk})« B 55 = B n 
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C ll “ WZ z! h y C 14 ~ p2 C /A 3’ C 22 C IV C 33 C 11 

C 41 = gezj(ph 3 ), c u = c u . = pz < /(pft 3 ), Cgg = C,, 

D 22 = 4 ^( 3 P^i) + D n = ^(p^) + P z < 2/ ( P^) 

D44 = P^(p^) + ^(spA*). D w = *(y - l)(4^ + xjfe* )/p 

£ 22= 4 ^( 3p ^) + M^(p^) + P^(p^). ^23= P^(3pA,A 2 ), £ 24 = py < /(3pA 1 A 3 ) 

B 32= ^( 3p M 2 ). £33 = £ 22> B 34= PV^( 3P V 3 )' E 42= P^( 3p Vl) 

E 43 = W(^3P*A). F 44 = E ll’ E 45 = (y - 1 )*(^ + ^ + ^)/P 

F 22 = = F 22’ F 44 = 4 ^( 3 P A 3 2 )' F 65=(y “ 1 )* V( P *5) 

The subscripts refer to the components of the rate-of-strain tensor. If the components are zero, they are not 
shown. 


7. Discussion of Results 


The method presented in Secs. 2 through 6 has been implemented into an operational code on a Univac 
82 system (equivalent to a Control Data Corporation (CDC) 6600) with a maximum storage of 263k words. The 
reason for selecting this aerobrake configuration for computational studies is that test data are available from 
the wind-tunnel and ballistic range, which are useful to check the method. The flow pattern observed in experi- 
ments indicates strong expansion over the shoulder, and the separating boundary layer on the back quickly 
merges into the near wake. A coupled solution of the shock layer in front of the nose and the flow surrounding 
the cylinder is deemed necessary to accurately predict the turning of the boundary layer and the base flow area. 

The code is designed to perform the flow simulation for a realistic spacecraft in the continuum regime be- 
tween 75 to 90 km and 9 to 12 kmps in altitude and speed, respectively. The major axis of the blunt frontal face 
has a magnitude between 6 and 10 m. The Reynolds number for this vehicle is in the range of 0.3 X 10^ to 0.5 
X 106 . The shock layer has a distinguishable viscous layer and a laminar near wake. For the present purpose 
of developing a computational method, a simple equation of state is used. The body is configured to have an 
ellipsoid with axis ratio 3:1, which is truncated by a flat disk and attached to a cylinder of radius equal to 1/3 
and length equal to the radius of the disk. A smooth contour is obtained in the mapped plane through four 
successive applications of Eq. (4). The hinge points, h y are located near the corners and edges of the aerobrake, 
shown in Fig. 2. The resultant mapped contour is closer to a pear shape than to a near circle; hence, the grid 
lines on top of the cylinder are still somewhat skewed. The origin of the physical polar coordinate system is 
located between the geometric center of the ellipsoid and the disk, so as to enable definition of the geometry by 
a single-valued function of 0. 

An initial grid of 10 x 36 x 7 (N x M x L) was first used to solve the Euler equations with as many as 300 
iterations. Then, a viscous grid of 28 x 36 x 7 was used to perform 400 more iterations. Despite the fact that an 
initial estimate for a 3D flowfield is always more difficult than for a 2D flowfield, the procedure given in Sec. 5 
is found to be quite satisfactory. At the end of 300 iterations, the Euler solution is generally converged in terms 
o(St and p w . Using the inviscid results to start the iteration for the viscous solution, changes of flow variables 
are mainly within the viscous layer and the near wake. The criteria of a converged viscous calculation is 
therefore based on the wall friction and heat-transfer coefficients, Cf and C/J, respectively. It has been noted 
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that these results tend to vary slightly from their mean values, and an accuracy level beyond the fourth-place 
digit can not be expected from the grid. More numerical experiments to investigate the effect of Cjv, damping 
parameters, and number of iterations have been reported in Ref. 9. The use of appended damping terms is 
critical to the success of the ADI computations. The procedures developed for blunt-body flow problems (Ref. 9) 
are found useful for the present finite-body problems. 

The computation time was approximately 6 hr central processing unit (CPU) for 100 iterations. The 
implicit procedure consumes 50% of the total time to set up the Jacobian matrices and to solve the block tri- 
diagonal equations. The mapping is estimated to use only 5% of the total time. The free-stream conditions are 
M = 34.8, Re = 8270, a = 18°, and Rn = 7 m. 

Figure 4 is a perspective view of the flowfield, wherein the upper portion of the pitch plane and the 
meridional plane at <}> = n/2 have the shock wave and Mach contours. A half of the body is shown in the 
background; some of the grid is removed for clarity. The shock location compares favorably with the schlieren 
photograph taken from the NASA Langley Research Center hypersonic wind tunnel at M = 10(1°) despite the 
different body geometries (Fig. 5). Another perspective view of the body is shown in Fig. 6, along with the 
trajectories of particles which move into the pitch planes on top of and behind the cylindrical body. Most of the 
particles follow the expected direction. The ring vortex and the crossflow downstream of the disk interact 
strongly and produce very intriguing flow patterns resembling a form of the trailing vortex behind a wingtip 
(Fig. 7). A closeup view of the grid is exhibited in Fig. 8. The grid is optimally distributed on the frontal face 
and packed toward the wall. A plot of Mach contours (Fig. 9) is used to study the shear-layer impingement on 
the lower portion of the cylinder. Sonic lines are plotted in dashed segments. Extensive recirculation areas are 
seen in Fig. 10 by means of velocity plot. The ring vortex starts at the sharp edge of the blunt nose and ends on 
the surface of the cylinder. The separation point of the boundary layer is further inward on the windward pitch 
plane as one would expect. The boundary-layer profile on the blunt face is also of strong interest, since it shares 
the same characteristics as observed in a low-speed plane stagnation flow. Lack of streamwise velocity 
gradient prevents development of the boundary layer as conventionally defined. Pressure and temperature 
contours are shown in Figs. 1 1 and 12, respectively, on the upper pitch plane. An expansion fan is observed 
between the body edge and the conic shock wave, and the near wake has a temperature range close to one-half 
the stagnation temperature. 

The viscous pressure distribution is given in Fig. 13. The flow incidence angle has shifted the distribu- 
tion to the negative x-axis. The base pressure is on the order of 10% the dynamic pressure. The pressure level 
on the upper disk is higher than that on the lower disk, because the flow speeds are different in those regions. 
Finally, the inviscid pressure distribution is given in Fig. 14 for comparative purpose. The edge effect is 
transmitted upstream and affects nearly 70% of the blunt face. 

8. Conclusion 


A method is proposed to reduce computational labor and storage requirement by taking advantage of the 
features in conformal mapping and the proven performance of an ADI factorization shock-fitting algorithm. In 
contrast to the methods based on the generalized coordinates, the number of metrics and the truncation error 
are minimized, the rate of convergence is substantially accelerated, and the level of accuracy is higher for the 
same number of grid points. The ADI technique is potentially competitive with the developing relaxation 
methods in the execution of large three-dimensional hypersonic problems on vector computers. In the contin- 
ual effort of certifying the present method, more quantitative comparisons with test data are being planned. 
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Fig. 1 Flowfleld schematics and physical 
coordinates. 
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Fig. 2 Sequence of coordinate transformation. 

a) Physical_space, X t Y,Z. b) Mapped 
space, X,Y,Z. c) Shear and cluster trans- 
formation, x,y f z. d) Computational space, 
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Fig. 9 Mach contours in the base flow region. 
(Dashed line = sonic line.) 


Fig. 10 Velocity plots in the base flow region. 


Fig. 11 Pressure contours in the upper pitch plane 
(normalized by p^). 


Fig. 12 Temperature contours in the upper pitch 
plane (normalized by T O0 ). 


X/XMflX 

X/XMFIX 

Fig. 13 Pressure distributions on the wall - NS 
solution. 

Fig.14 Pressure distributions on the wall - Euler 
solution. 
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